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ABSTRACT 

We present nonlinear three-dimensional simulations of the stably-stratified portion of the solar 
tachocline in which the rotational shear is maintained by mechanical forcing. When a broad toroidal 
field profile is specified as an initial condition, a clam-shell instability ensues which is similar to the 
freely-evolving cases studied previously. After the initial nonlinear saturation, the residual mean fields 
are apparently too weak to sustain the instability indefinitely. However, when a mean poloidal field is 
imposed in addition to the rotational shear, a statistically-steady state is achieved in which the clam- 
shell instability is operating continually. This state is characterized by a quasi-periodic exchange of 
energy between the mean toroidal field and the instability mode with a longitudinal wavenumber 
m = 1. This quasi-periodic behavior has a timescale of several years and may have implications for 
tachocline dynamics and field emergence patterns throughout the solar activity cycle. 
Subject headings: Sun: tachocline, MHD, Sumrotation, Suminterior 



1. INTRODUCTION 

It is now well established (at least from a theoreti- 
cal perspective) that axisymmetric rings of toroidal field 
in the solar tachocline are susceptible to global insta- 
bilities induced by the latitudinal differential rotation. 
This was first demon st rated f or 2-D spherical s u rfaces 
dGilman fc Foxl [l997l fl99l IDikpati fe Gilmanl [19991: 
Gilman fe Dikpatill200q ICallv et al.ll2003t iDikpati et al l 
2004al) and was later extended to 3-D spherical shells un- 
der the shallow-water app roximation (|Gilman fc D ikpati 
l2002t Dikpati et al. 200 3) and the thin-shell approxima - 
tion ((CallyJ 2003; Miesch et aLll2007t iGilman et al.ll2007f L 
For strong fields (such that the magnetic energy is com- 
parable to the kinetic energy contained in the differen- 
tial rotation) , the most unstable modes have longitudinal 
wavenumber m = 1, which corresponds to a tipping of 
the ring such that its central axis becomes misaligned 
with the rotation axis. 

This class of global magneto-shear instabilities occurs 
both for concentrated bands of toroidal field and for 
broad field profiles. The preferred mode for broad field 
profiles which are antisymmetric about the equator is 
the clam-shell instability whereby toroidal rings in the 
northern and southern hemispheres tip out of phase, re- 
connecting at the equator on o ne side of th e sphere and 
openi ng up on th e othe r side (fCalTvl I200H: ICallv et all 
120031 ; IMiesch et aTl 2007'). In the absence of external forc- 
ing, the instability proceeds until loops of field become 
perpendicular to the equatorial plane and the latitudinal 
shear is nearly eliminated. 

However, the solar tachocline is not an isolated system. 
Rotational shear is continually maintained by stresses 
from the overlying convection zone and the underlying 
radiative interior. Meanwhile, magnetic flux is continu- 
ally being replenished from above by penetrative convec- 
tion and meridional circulation. If clam-shell instabilities 
are indeed occurring in the tachocline, their environment 
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is likely much more dynamic than the freely-evolving 
scenarios considered in previous nonlinear simu lations 
(|Callvll2001l : ICallv et al.ll200l IMiesch et al.ll2007h . 

In this letter we present the first nonlinear simula- 
tions of global magneto-shear instabilities in the so- 
lar tachocline in which the instabilities are continually 
maintained against nonlinear saturation and dissipation 
through the use of external forcing. We use the same non- 
linear , 3-D thin-shell model employed by IMiesch et alj 
(|2007l hereafter MGD07) to study the freely-evolving 
case. The linear st ability of this t hin-sh ell system has 
been investigated by I Gilman et a l. (2007). The reader is 
referred to these two papers for a much more compre- 
hensive discussion of the unstable modes and their non- 
linear saturation under freely-evolving conditions. These 
papers also contain a more comprehensive discussion of 
related work and a more detailed description of the thin- 
shell model and the numerical algorithm. 

The first question we address in this letter is whether 
the instabilities can generate poloidal fields efficiently 
enough to act as a dynamo when only mechanical forc- 
ing of the rotational shear is applied. Our simulations 
suggest that this is unlikely although the results are in- 
conclusive. However, statistically steady states are found 
when the magnetic energy is maintained via a poloidal 
forcing term in the mean induction equation. These 
states exhibit a quasi-periodic behavoir as mean fields 
continually transfer energy to the non-axisymmetric in- 
stability modes and are then re-established by the me- 
chanical and magnetic forcing. 

2. THE THIN-SHELL MODEL 

Our numerical model is based on the thin-shell approx- 
imatio n which is described in detail bv lMiesch fc Gilmanl 
(2004) and MGD07. The aspect ratio S = D/R is as- 
sumed to be much less than unity, where D and R are 
the width and radial location of the computational do- 
main. Only terms of lowest order in S are retained in 
the equations of motion and the layer is assumed to be 
stably-stratified, cor responding to the lowe r portion of 
the solar tachocline ([Thompson et al.ll2003h . The equa- 
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tions are expressed in a rotating, spherical polar coordi- 
nate system denned by the colatitude 9, longitude <j> and 
height z with corresponding unit vectors 6, 4>, and z. 
The velocity and magnetic field components are defined 
as v = u<f>+v9+Swz and B = a(f)+b0+6cz. The system 
is made nondimcnsional through the use of horizontal 
and vertical length scales R and D, a velocity scale based 
on the equatorial rotation rate, and the background den- 
sity and entropy gradient. We neglect pressure variations 
relative to entropy variations in the linearized equation 
of state which corresponds to the Boussinesq limit dis- 
cussed by MGD07 [8p = in their equation (6)]. Apart 
from the parameters which specify the forcing, dissipa- 
tion, and initial conditions, this leaves one free param- 
eter: the Froude number F r , which is a nondimcnsional 
measure of the stratification. For the simulations dis- 
cussed here we take F r = 0.1 as appropriate for the lower 
tachocline (MGD07). 

The subgrid-scale (SGS) model consists of a fourth- 
order hypcrdiffusion in the horizontal dimensions with 
an effective Reynolds number of Rh and a Smagorinsky 
formulation for the vertical diffusion involving a verti- 
cal Reynolds number R v (MGD07). The amplitude and 
form of the SGS magnetic and thermal diffusion are as- 
sumed to be the same as the viscous diffusion, implying 
magnetic and thermal Prandtl numbers of unity. The 
upper and lower boundaries of the layer are assumed to 
be impenetrable, stress-free, perfectly conducting, and 
isentropic. 

Mechanical forcing is imposed by adding a volumetric 
source term to the momentum equation of the form 

F u {6,z)=t- 1 {u -{u))4> (1) 

where u is the zonal velocity, r is a characteristic 
timescale for the establishment of the shear and 

uq(9, z) = — (1 — cos (nz)) (so — s 2 cos 2 9) sin# . (2) 

Angular brackets <> in equation |T]) denote an average 
over longitude. In equation sq is the nondimcnsional 
equatorial rotation rate relative to the rotating reference 
frame and s 2 is the fractional angular velocity difference 
between the equator and poles. Here we take so = 0.044 
and s 2 = 0.18 (MGD07). The vertical profile of u is 
such that the latitudinal shear is maximum at the top of 
the layer (z = 1) and vanishes at the bottom (z = 0) in 
analogy with the solar tachocline. 

In selected cases we also impose a dipolar field by 
means of a similar forcing term added to the induction 
equation: 

F p (0,z)=t- 1 (Bo-(B p )) (3) 

where B p = bO + Scz and 

B Q (0, z) = VxVx [C o cos0sin(7rz)£] . (4) 

This field is intented to represent generated by dynamo 
processes in the convection zone and pumped downward. 
Imposing Bq together with uq is an indirect way of im- 
posing a mean toroidal field; the latitudinal shear in the 
upper portion of the layer stretches and amplifies the 
poloidal field through what is known as the fi-effect, re- 
sulting in a mean toroidal field (a) which is antisymmet- 
ric about the equator. The peak amplitude of (a) will in 
general depend on the poloidal field strength Co, on the 
SGS diffusion, and to some extent on F r . 
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Fig. 1. — The volume-integrated magnetic energy contained in 
the mean toroidal field (TFME, solid line), the mean poloidal field 
(PFME, dash-dotted line) and the non-axisymmetric field (NAME, 
dotted line) are shown as a function of time in a simulation with 
no magnetic forcing (R^ = 10 7 , R v = 10 4 , r = 0.1). All curves 
are normalized relative to the integrated kinetic energy contained 
in the imposed shear profile uq (indicated by the horizontal dashed 
line). 
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Fig. 2. — As in Figure IT1 for a case with both mechanical and 
magnetic forcing as described in the text (R^ = 10 5 , R v = 10 4 , 
t = 0.1). Note that the vertical axis is linear here rather than 
logarithmic. Vertical dashed lines indicate the time instances illus- 
trated in Figure [3] 

3. SUSTAINED MAGNETO-SHEAR INSTABILITIES 

Our simulations are initiated as described in MGD07, 
with an equilibrium state defined by the zonal velocity 
Uq and a broad toroidal field profile of the form 

a(9) = a cos sin 9 . (5) 

The initial pressure and temperature are chosen such 
that the initial state is in magneto-hydrostatic and 
magneto-geostrophic balance. The parameter a is set 
to unity for the simulations presented here, correspond- 
ing to a peak dimensional field strength of about 40 kG, 
toward the high end of the range expected to exist in 
the tachocline (MGD07). Global magneto-shear insta- 
bilities also occur for weaker fields but they take longer 
to develop. The initial equilibrium state is perturbed by 
adding a random, small-scale velocity field. 

Unlike the simulations presented in MGD07, we main- 
tain the differential rotation through the forcing term 
defined in equation ([1]). An example of the subsequent 
evolution in the absence of magnetic forcing (F p = 0) is 
shown in Figure [T] The spatial resolution used for this 
case is N e , iV , N z = 128, 256, 210. 



Magneto-Shear Instabilities in the Tachocline 



3 



i 



1*1 IPO 




Fig. 3. — The horizontal structure of the magnetic potential J near the top of the layer (z = 0.98) is illustrated at three different times 
for the same simulation as shown in Figure [2] The color table varies with increasing J as blue-green-yellow-red. Black countours trace 
horizontal magnetic field lines as noted in the text, although the sense of the field is opposite in the northern and southern hemispheres. 
Each horizontal surface is displayed as a Molleweide projection in which lines of constant latitude appear horizontal. 



Figure [T] shows the components of the magnetic energy 
integrated over the volume of the shell. One time unit 
in the nondimensional system corresponds to about four 
days (MGD07). At early times the mean toroidal field 
(TFME) dominates but the non-axisymmetric compo- 
nent (NAME) grows rapidly as the clam-shell instability 
develops. Most of the NAME is in the m = 1 mode which 
dominates the total magnetic energy after the instability 
saturates at t ~ 400. 

The evolution shown in Figure [1] is very similar to the 
unforced cases discussed at length in MGD07. However, 
in the absence of mechanical forcing, the saturation of 
the clam-shell instability induces a global redistribution 
of angular momentum which reverses the sense of the 
differential rotation (MGD07). The forcing suppresses 
this, leaving the differential rotation profile unchanged 
after saturation. 

After saturation, the magnetic energy in the mean 
toroidal field decreases as in the unforced case, despite 
the persistent rotational shear. The TFME oscillates 
with a period of about 15 time units, corresponding to 
about two months. At later times, the amplitude of the 
oscillation decreases and the period increases slightly, to 
about 18 time units. This oscillation appears to be in- 
duced by a standing Alfven wave excited by the initial 
saturation of the instability. However, low-wavenumber 
Rossby waves are also excited and have a comparable 
period. 

In shallow-water systems which have a deformable 
upper boundary, global magneto-shear instabilities can 
possess significant kinetic helicity, suggesting they may 
serve to generate poloidal field from toroidal field 
and thus drive a self-sustain ed dynamo contained en- 
tirely within the tachocline (iDikpati fc Gilmanl l2001at 
iGilman fc Dikpatill2002t IDikpati et al.ll2003D . 

The simulation shown in Figure[T]does not extend more 
than a diffusive timescale (~ 10 5 time units) so it is un- 
certain whether or not it may be classified as a dynamo. 
The energy in the mean fields may be levelling off be- 
yond t ~ 800 but it is unclear whether this state will 
persist indefinitely. A more diffusive analogue of this 
case (Rh = 10 5 ) was certainly not a dynamo since the 
total magnetic energy decayed steadily after the initial 
saturation. 

In any case, it appears either that the clam-shell in- 
stability is no longer operating beyond t ~ 500 or that it 
is operating on a much longer time scale. The magnetic 
energy remains dominated by the m = 1 component for 
at least five years in dimensional time units. Over this 
time period the energy in the mean toroidal field remains 
larger than the equipartition energy associated with the 



imposed differential rotation but the energy in the mean 
poloidal field is about three orders of magnitude smaller. 
These fields appear to be a remnant of the initial satura- 
tion of the instability at t ~ 400 as opposed to dynamo- 
generated fields induced by ongoing instabilities. 

The situation changes dramatically if a mean poloidal 
field is imposed through the forcing term expressed in 
equation (J3]) . Figure [2] illustrates the evolution of the 
magnetic energy components in a simulation with both 
mechanical and magnetic forcing (Ng = 64, N$ = 128, 
N z = 210). The time span shown covers a period beyond 
the initial saturation of the instability, after a statisti- 
cally steady state has been reached. 

The time evolution shown in Figure [21 reflects a quasi- 
periodic exchange of energy between the mean toroidal 
field and the non-axisymetric field components, the latter 
dominated by the m = 1 mode. The first two oscillations 
shown each span about 340 time units which corresponds 
to 3.7 years. However, the TFME drops lower in the sub- 
sequent cycle, reaching a minimum at t ~ 1950, leading 
to a longer cycle of about 400 time units (4.4 years). 
The following cycle is then much shorter, lasting only 
about 230 time units (2.5 years). The shape of each cy- 
cle is asymmetric, with a relatively slow rise in TFME 
followed by a sharper drop as the clam-shell instability 
sets in. 

The non-axisymmetric magnetic energy exhibits quasi- 
periodic cycles similar to the mean toroidal field but 
phase-shifted such that the maxima in NAME occur as 
the TFME is decreasing. Again, this reflects the repeated 
development of the clam-shell instability which transfers 
energy from the mean toroidal field to the m = 1 com- 
ponents. After the instability saturates, the mechanical 
and magnetic forcing re-establish the mean fields and the 
next cycle proceeds. 

To gain further insight into the dynamics, it is instruc- 
tive to express the magnetic field in terms of scalar mag- 
netic potentials J and C defined such that 



B = Vx (Jz) + Vx Vx (Cz) 



(6) 



In the simulations reported here, as in the freely evolving 
cases reported in MGD07, the magnetic field remains 
predominantly horizontal and the first term on the right- 
hand-side of equation |6]) dominates. Thus, to a good 
approximation, contours of J trace the horizontal field 
lines as illustrated in Figure [3] [V X (Jz) = zxVJ]. 

The changing patterns shown in Figure [3] illustrate the 
competing effects of the forcing and the instabilities. At 
t = 1800 (Fig. [3)2), the mean toroidal field dominates the 
magnetic energy, although non-axisymmetric structure 
is evident. By t = 1950 the clam-shell instability has 
transferred much of this energy to the m = 1 mode and 
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horizontal field lines are oriented more north-south (Fig. 
E3>). The imposed shear then operates on these fields as 
well as the imposed poloidal field to rebuild the mean 
toroidal field which again dominates by t = 2150 (Fig. 

m- 

The vertical structure of the flow is similar to analo- 
gous freely-evolving cases with vertical shear discussed 
in MGD07. The velocity and magnetic fluctuations re- 
main predominantly horizontal and the instability pro- 
ceeds most vigorously near the top of the layer where the 
latitudinal shear is strongest. 

For the simulation shown in Figures [2] and [3l the am- 
plitude of the imposed poloidal field, Co, is such that 
the integrated magnetic energy PFME is about twice 
the equipartion value Uq/2 (integrated over the volume). 
This ratio is somewhat smaller near the top of the layer 
where the latitudinal shear peaks but nevertheless it is 
probably unrealistically large for the solar tachocline. 
Weaker imposed fields could potentially produce mean 
toroidal fields of comparable strength but this can only 
be achieved in a simulation if the SGS diffusion is suf- 
ficiently low. Indeed, an analogous simulation with a 
weaker imposed field (PFME/DRKE = 0.02) and the 
same diffusion coefficients (Rh — 10 5 , R v — 10 4 ) pro- 
duced weaker mean toroidal fields (TFME/DRKE - 3). 
This system also exhibits sustained clam-shell instabili- 
ties with quasi-periodic behavior but the phase relation- 
ship between the TFME and NAME is not as well de- 
fined as in Figure [2] Some quasi-periodic oscillations are 
present on timescales of several years but longer-term 
trends are also evident, comparable to or longer than 
the duration of the simulation (~ 10 years). Such longer- 
term evolution is to be expected since the growth rate of 
the clam-shell instability decreases with decreasing field 
strength. 

Achieving substantially lower SGS diffusion would re- 
quire higher resolution which is a challenge because of 
the long integration times necessary to capture multiple 
cycles. The relatively short forcing timescale used for the 
simulations shown in Figures [T][3](t = 0.1, corresponding 
to about 10 hours) is also in some sense required by the 
setup of our numerical experiments. In a simulation sim- 
ilar to that shown in Figures [2][3] but with r = 2 (8 days), 
the mechanical forcing is insufficient to overcome the 
magnetic tension associated with the imposed poloidal 
field. As a result, the DRKE is only 20% of the target 
value associated with uq and the clam-shell instability is 
suppressed. The field does exhibit m = 1 structure but 
TFME dominates the magnetic energy and the evolution 
is quasi-steady, with a slow retrograde propagation and 
no episodic opening up of the clam-shell pattern. 

The solar tachocline is much more complex than our 
highly simplified model. Forcing timescales are longer 
but the diffusion is lower, so strong toroidal fields could 
be produced with lower poloidal field strengths and the 
differential rotation could be maintained with weaker 
mechanical forcing. It is therefore difficult to establish 
from the simulations alone whether the solar tachocline 
is in a regime which exhibits quasi-periodic behavoir as 
in Figure [21 However, helioseismic inversions indicate 
that the rotational shear in the tachocline is continually 
maintained and linear analysis suggests that even weak 
toroidal fields are unstable in the presence of such shear 



so it is likely that the tachocline is indeed continually 
undergoing global magneto-shear instabilities. 

Latitudinal shear in the convective envelope is main- 
tained by Reynolds stresses and meridional circulations 
on timescales of order months while the nearly uniform 
rotation of the radiative inter ior may be established on 
much longer time scales (e.g. iGough fe Mclntvrd 119981 : 
iTalon et al]|2002HMieschll2005l ). Magnetic flux is contin- 
ually supplied to the tachocline by topological pu mping 
from penetrative convection fe.g. lTobias e t al. 200rj). Al- 
though much of this flux is disordered, mean fields may 
be generated by nonlinear self-organization processes 
such as inverse cascades or by rotational phase mixing 
and t urbulent reconnection (|Spruitlll999l : lBrowning et al.l 

2006) . In addition, global meridional circulations may 
transport axisymmetric poloidal flux to the tachocline 
where mean toroid al fields would then be generated 
by ro tational shear (|Dikpati et alj l2004bt ICharbonneaul 
I2005T) . 

In this letter we have demonstrated that the clam- 
shell instability can operate continually when rotational 
shear and magnetic fields are continually replenished. 
The temporal evolution is quasi-periodic and as mean 
fields alternately build up and destabilize. This has 
the character of a critical phenomenon but it is not 
self-organized criticality in the technical sense because 
there is a characteristic time and spatial scale associated 
with the grow th rate and wavenumber of the instability 
(|Jensen|[l998f) . It is likely that other instability modes 
may similarly be maintained continually by external forc- 
ing. Notable among these is the m = 1 tipping insta- 
bility which exists for the type of concentrated toroidal 
bands thought to give rise to photospheric active regions 
(IDikpati fc Gilman|[l999tlCallv et atll2003t [Miesch et all 

2007) . 

If global magneto-shear instabilities are indeed occur- 
ring in the solar tachocline they would have wide-ranging 
implications for tachocline dynamics and for the coupling 
between the convective envelope and the radiative inte- 
rior. Angular momentum transport induced by these in- 
stabilities may influence the differential rotation profile 
in the convection zo ne and the longer-term rotational 
evolution of the Sun (ICharbonneau fc MacGregolll993t 
iGilmanl l2000t ITalon et al.l 120021 ). Chemical transport 
across the tachocline also has implications for solar evolu- 
tion, helioseismic structura l inversions, and photospheric 
abundance measure ments (jChristensen-Dalsgaardl 120021 : 
iPinsonneaultl [1997). Tachocline instabilities may also 
play a role in the solar dynamo. Although the simu- 
lations reported here suggest that clam-shell instabili- 
ties may not be capable of sustaining a dynamo local- 
ized entirely within the tachocline, they may still play a 
role in parity se lection, field propagation, and flux emer- 
gence patterns ([Gi lman 200 (1: IDikpati fc GilmanlfeOOlbl ; 
iNorton fc Gilmanll2005l : IDikpati fc Gilmanll2005fl . 
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